Effect of particle statistics in strongly correlated two-dimensional Hubbard models 
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We study the onset of particle statistics effects as the temperature is lowered in strongly correlated 
two-dimensional Hubbard models. We utilize numerical linked-cluster expansions and focus on the 
properties of interacting lattice fermions and two-component hard-core bosons. In the weak-coupling 
regime, where the ground state of the bosonic system is a superfluid, the thermodynamic properties 
of the two systems at half filling exhibit very large differences even at high temperatures. In the 
strong-coupling regime, where the low-temperature behavior is governed by a Mott insulator for 
either particle statistics, the agreement between the thermodynamic properties of both systems 
extends to regions where the antiferromagnetic (iso)spin correlations are exponentially large. We 
analyze how particle statistics affects adiabatic cooling in those systems. 
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I. INTRODUCTION 

The Fermi-Hubbard model has been the de facto play- 
ground for exploring the properties of high-temperature 
superconductors for more than two decades Yet, still 
no analytical solution exists in more than one dimen- 
sion, and state-of-the-art numerical calculations prove 
very difficult in regimes where the average number of 
fermions per site is different from 1. Unveiling the prop- 
erties of the model in the latter regime, and addressing 
whether it supports superconductivity, may be crucial in 
understanding high-temperature superconductivity Q. 

More than a decade ago, it was proposed that one could 
"solve" strongly-interacting quantum lattice models by 
"simulating" them using ultracold atoms in optical lat- 
tices 3. More recently, the Mott insulator in the Bose-Q 
and the Fermi-Hubbard || models were realized in these 
experiments. However, current accessible temperatures 
for fermions are still higher than one needs to observe 
even the relatively high temperature antiferromagnetic 
Neel transition in three dimensions. 

These experiments are done using atoms with inter- 
nal degrees of freedom, which can be selected to emulate 
not only two-specie fermions (Fermi-Hubbard model) and 
single-specie bosons (Bose-Hubbard model) , but also par- 
ticles with exotic statistics and/or pseudo-spins [oHlOf. 
For instance, it has been shown that experiments with 
two-component (spin-1/2) bosons can lead to the real- 
ization of quantum spin models with tunable parame- 
ters [Tl| - [T5j . In the specific case where the intra-specie 
onsite repulsion is infinite, multiple occupancy of a sin- 
gle specie per site is forbidden. This case realizes an 
effective two-component hard-core boson (2HCB) model, 
which could be thought of as the bosonic equivalent of the 
Fermi-Hubbard model [IH EH, [l?} • Similar to fermions, 
in the strong-coupling regime (large inter-species inter- 
actions), the low-energy properties of this model can be 
described by a t-J model pH [l9l| . 

Despite the outward similarities between the Fermi- 
Hubbard model and the 2HCB-Hubbard model (in one 



dimension they share identical thermodynamic proper- 
ties) in two dimensions, particle statistics plays a funda- 
mental role in the properties of the system as the tem- 
perature is lowered, and in the selection of the ground 
state. At half filling, the ground state of fermions in two 
spatial dimensions is a Mott insulator with a long-range 
Neel order for any value of the interaction strength, U. 
For 2HCBs, there is a quantum critical point at interac- 
tion U c /t ~ 11, where t is the hopping amplitude, which 
separates a phase with two miscible strongly-interacting 
superfluids (2SF) at small inter-specie interactions from 
a Mott insulator super-counter-fluid (SCF) [l]| phase in 
the strong-coupling regime [ljj]. The latter state corre- 
sponds to a superfluid of pairs of bosons from one specie 
and holes of the other specie, and can also be interpreted 
as a long-range YY-ferromagnet in the iso-spin language. 
It is expected that such big contrasts in the nature of 
the ground states result in significant differences in the 
thermodynamic properties as well. Hence, for validating 
experiments with ultracold gases, which are performed at 
finite temperature, it is important to have access to exact 
quantitative results for the finite-temperature properties 
of the two systems. 

While there have been numerous finite-temperature nu- 
merical studies of fermions in optical lattices 0, [2l| . 
the same is not true for 2HCBs, in particular in two 
dimensions. In one dimension, calculations were done 
introducing a generalized Jordan- Wigner transforma- 
The three-dimensional model with attractive 
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tion 

interaction was studied to describe the supersolid state 
of 4 He [23|, [24|. Recently, dynamical mean-field the- 
ory results have been reported for the two-component 
(soft-core) Bose-Hubbard model in two dimensions with 
an average of 1/2 particle of each specie per site and 
large intra-specie interactions [TBI . [25j . It was found 
in those studies that, upon heating from zero tempera- 
ture, the system quickly enters an unordered Mott insu- 
lator. These findings have been complemented by finite- 
temperature quantum Monte Carlo (QMC) simulations 
of magnetic phases [13, as well as a field-theoretical 
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treatment [27j in the hard-core limit. 



II. THE MODEL AND NLCES 



Here, we utilize numerical linked cluster expansions 
(NLCEs) to provide a comparative analysis between 
finite-temperature properties, such as the equation of 
state, entropy, specific heat, double occupancy, and spin 
correlations, of fermions and 2HCBs. We are particu- 
larly interested in identifying at what temperatures par- 
ticle statistics become important in different interaction 
regimes, as well as what kind of qualitatively different be- 
havior is produced by the statistics of the particles below 
those temperatures. We also discuss the implications of 
having bosonic vs fermionic statistics for adiabatic cool- 
ing protocols and for detecting short-range spin correla- 
tions in ultracold atoms experiments. We note that the 
lowest temperatures that are accessible with NLCEs are 
typically higher than the crossover temperatures to the 
Xy-fcrromagnet phase [15) . Therefore, the behavior of 
the strongly-correlated systems that we study at half fill- 
ing is that of a Mott insulator with large spin correlations 
for either particle statistics. 



We consider the two-dimensional (2D) Hubbard Hamil- 
tonian on the square lattice: 



ala ja + K.c.) + Uj2' 



(1) 



where a| {(Xi a ) creates (annihilates) a particle with spin 
a (for simplicity, we use spin instead of isospin for bosons 
too) on site i, and — a\ a a,i a is the number operator. 
(..) denotes nearest neighbors (NN), and U (> 0) is the 
strength of the onsite repulsion, t = 1 (h = 1 and ks = 
1) sets the unit of energy throughout this paper. We 
consider two different particle statistics: fermions (a^ = 
fia) and 2HCBs (ai a = bi a ). The 2HCB operators satisfy 
the following commutation relations and constraints: 



m = hi=o. 



(2) 



We show that, in the weak-coupling regime (U < U c ), 
where the bosonic system has a superfluid ground state, 
the disagreement between properties of fermions and 
2HCBs is apparent at relatively high temperatures (of 
the order of the hopping). In contrast, by increasing 
U in the strong-coupling regime, the agreement between 
the thermodynamic properties of the two systems ex- 
tends to lower temperatures. For 2HCBs, although the z- 
antifcrromagnetic Neel ground state is known to win over 
the Xy-ferromagnet only in cases where the two species 
have different hopping amplitudes, the two phases are 
very close in energy for equal hopping amplitudes (the 
case considered here) as only terms beyond the second 
order perturbation in interaction determine the differ- 
ence [ll|, [llj]. Consistent with this picture, we find that 
at strong interactions, short-range z-antiferromagnetic 
correlations are large in the low-temperature Mott re- 
gion and very close for the two particle statistics. Given 
that there are more efficient cooling techniques for bosons 
than for fermions, one could envision probing finitc- 
tcmpcrature fermionic correlations using strongly corre- 
lated bosonic systems. Here, we present evidence that, 
through an adiabatic cooling mechanism that takes place 
in the bosonic system, the region with exponentially large 
antifcrromagnetic (AF) correlations in two dimensions is 
more easily accessible with 2HCBs than with fermions. 



The exposition is organized as follows: In Sec. [TJl we 
introduce the model and discuss the NLCEs used in this 
study. We present the results for the thermodynamic 
properties of the models as well as their implications for 
the optical lattice experiment in Sec. IIII1 Our findings 
arc summarized in Sec. IIVI 



Numerical linked-cluster expansions 

We solve the Hamiltonian ([1]) using the NLCEs intro- 
duced in Ref. [28) . In NLCEs, an extended property of 
the lattice model per site in the thermodynamic limit, 
P, is expanded in terms of contributions from all of the 
clusters, up to a certain size, that can be embedded in 
the lattice: 



P = ^L(cK(c), 



(3) 



where c represents the clusters. The contribution of each 
cluster with a particular topology is proportional to the 
number of ways it can be embedded in the lattice per 
site, L(c), and its weight for the property of interest, 
Wp{c). The weights are computed based on the inclusion- 
exclusion principle and given the property for each clus- 
ter, 'P(c), which is calculated using exact diagonaliza- 
tion 1281: 



w p (c) = V{c) - ^w p { 



(4) 



Here, we carry out the calculations up to the ninth order 
in the site expansion (maximum cluster size of nine sites). 

NLCEs do not suffer from statistical or systematic er- 
rors, such as finite size effects and are not restricted to 
small or intermediate interaction strengths. For this rea- 
son, they are complementary to more commonly used 
methods, such as QMC simulations and dynamical mean- 
field theory, especially in the strong-coupling regime 
(U^>t) where computations in the latter approaches be- 
come more challenging. However, NLCE results are use- 
ful only in the temperature region in which the series con- 
verge, which has been shown to extend beyond the region 
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accessible within high-temperature expansions [28|, l29j . 
The convergence of NLCEs can be accelerated by means 
of numerical resummations algorithms pjjj ]. Here, we use 
Euler and Wynn methods with different parameters and 
plot the resulting last two orders (or only the last order 
when the two are indistinguishable in the figures). The 
results from all of those algorithms agree with each other 
in the regions shown within the small fluctuations seen 
in some cases at the lowest temperatures. These fluctu- 
ations, which occur below the convergence temperature 
of NLCE direct sums, arise from numerical instabilities 
in the resummations routines. 

We work in the grand canonical ensemble [29j], and so, 
the exact diagonalization for every finite cluster in the se- 
ries is performed in all particle and spin sectors. For each 
U, the partition function and all other observables are 
calculated in a dense grid of chemical potentials (/i) and 
temperatures (T). This allows us to study their behavior 
at constant density n = (n-|- + %) = 2(h a ) [29l ]. where (..) 
denotes the expectation value. Since only NN hopping is 
considered, properties of the particle-doped system can 
be expressed in terms of those for the hole-doped system. 
Hence, in most cases away from half filling, we show re- 
sults only for the hole-doped system (n < 1). 



III. RESULTS 




Entropy and specific heat 

Generally, when using QMC-based methods, the spe- 
cific heat (C v ) or the entropy (S) calculations involve nu- 
merical derivatives and/or integration by parts [30L [3l|. 
which can introduce systematic errors. Within NLCEs, 
these two quantities are computed directly from their def- 
initions: 



S = ln(Z) + 



(H) - »(n) 
T 



(5) 



where Z is the partition function, and 



' d(H) \ (d(H)\ (d(H) 



df 



(6) 

Since we work in the grand canonical ensemble, where 
the chemical potential and not the density is the control 
parameter, we have written this expression in a more suit- 
able form for numerical evaluation. After straightforward 
mathematical derivations, using Maxwell equations, one 
can obtain the following closed form for the specific heat 
in terms of expectation values that can be computed di- 
rectly in NLCEs: 



Cy 



1 

T 2 



(AH 2 ) 



((Hn)-(H)(n) 
(Af^j 



(7) 



FIG. 1. (Color online) (a) Specific heat (C v ) of the Fermi- 
Hubbard model at half filling vs temperature for values of 
the onsite interaction ranging from U — 2 to three times the 
bandwidth (U = 24). (b) Difference between the specific heat 
of the Fermi-Hubbard model and the 2HCB-Hubbard model 
at half filling vs temperature, for the same interactions as in 
(a). The left inset in (b) shows C v for the latter model vs T 
for U = 2, 4, 6 and 8; and the right inset shows the difference 
between C v of the Fermi-Hubbard model and 2HCB-Hubbard 
model at half filling vs U at three different temperatures. We 
have used Euler sums for the last 6 terms. Thick (color) lines 
are the results of the sums up to the 9th order and thin (black) 
lines up to the 8th order. 



where (AH 2 ) = (H 2 ) - (H) 2 . and similarly (An 2 ) = 
(n 2 ) - (h) 2 . 

We begin our study of the dependence of the observ- 
ables on particle statistics as the temperature in the sys- 
tem is changed by showing, in Fig. [TJ the specific heat 
for fermions and 2HCBs. Figure [Ha) shows the specific 
heat of the Fermi-Hubbard model from the weak-coupling 
regime to the strong-coupling regime with interactions 
up to three times the bandwidth (U = 24). The trend 
in the deviation of the specific heat for 2HCBs from that 
of fermions for different interaction strengths is depicted 
in Fig. [ljb) , where we show the difference between the 
two. It is clear that for U < 8, the difference is signif- 
icant at temperatures greater than 1. For example, as 
shown in the left inset in Fig. HJb) , the specific heat of 
2HCBs for U — 2 is roughly twice as large as that of 
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FIG. 2. (Color online) (a) Temperature and (b) value of the 
high-temperature peak of the specific heat (C^) as functions 
of U. By decreasing U in the strong-coupling regime, Cf 
for the bosonic case deviates from the fermionic one around 
U — 12. For U < 8, the two peaks in the specific heat of 
2HCBs merge, and the high-temperature peak is not well- 
defined anymore, (c) The location of the high-temperature 
crossing point between the specific heat curves for consecutive 
values of U. (d) The value of C„ at that crossing point. 



fermions around T = 1 . Although the exact trend of the 
former at lower temperatures cannot be resolved within 
our method, such large high-temperature values in com- 
parison to those for the fermionic case are suggestive of 
the lack of a second peak at lower temperatures in the 
2HCB model. This is also supported by a fast drop in 
the entropy (not shown). We have found that for (7 = 2, 
at T ~ 0.8, already ~ 70% of the infinite-temperature 
entropy has been quenched. For 2HCBs, a double-peak 
structure in the specific heat becomes apparent for U > 7. 
At the same time, the minimum temperature at which 
NLCEs converge, which for smaller U is generally higher 
for 2HCBs in comparison to fermions, extends to roughly 
the location of the low-temperature peak. 

For fermions, our NLCE calculations resolve the 
double-peak structure of the specific heat in the ther- 
modynamic limit for U > 6. Since QMC simulations 
can access lower temperatures for U < 6, previous QMC 
studies of this model have established the existence of 
the low-T peak for any finite value of the interaction 
strength j32[, which signifies the crossover to the phase 
with exponentially large AF correlations. In the strong- 
coupling regime, the specific heat results for fermions and 
2HCBs are very close to each other at high temperatures 
(T > 1). In this regime, the difference between the C v of 
the two particle statistics decreases systematically for all 
accessible temperatures as U is increased. As shown in 
the right inset in Fig.QJb), the reduction of this difference 
is nearly exponential for large values of U. 

The high-T peak, which is associated with the frccz- 




(c) u=16 



10 0.1 



10 



FIG. 3. (Color online) Comparison of the entropy vs tem- 
perature for n = 0.85 (thick lines) and n = 1.00 (thin lines) 
between the Fermi-Hubbard and the 2HCB-Hubbard models 
for (a) U = 4, (b) U = 8, and (c) U = 16. 



ing of charge degrees of freedom and moment formation, 
moves to higher T as U increases. In the atomic limit 
(t — > 0), the location of the peak, T p , is determined by 
atanh(a) = 1, where a = C7/16T P , which can be ap- 
proximated by T p = C//4.8. For large U, and regardless 
of the particle statistics, we find a very good agreement 
with the atomic limit prediction for T p [see Fig. (2{a)]. 
Note that there is no well-defined high-T peak in the spe- 
cific heat of 2HCBs for U < 7. Unlike its position, the 
value of this peak does not change monotonically with 
increasing U and, as seen in Fig. HJb), has a minimum 
around U = 7 for fermions and around U — 9 for 2HCBs. 
However, the relative change in the studied range of inter- 
actions is only about 20%. For large values of U, the area 
under the high-T peak of C v /T also approaches ln(2), 
consistent with results in the atomic limit. 

The low-temperature peak, which is associated with 
ordering of the moments, is expected to move to lower 
temperatures by increasing U. This is because, for large 
values of C7(> 12), the low-T system is essentially de- 
scribed by the antiferromagnetic Heisenberg model with 
a characteristic energy scale of J oc t 2 /U [21(. Therefore, 
the position of the low-T peak, which we do not report 
here, is expected to be inversely proportional to U . 

An interesting feature discussed in the past for cor- 
related fermionic systems is the near universal higlr-T 
crossing of the specific heat curves for different values 
of the interaction [32T - l34j ]. Due to its accuracy, NLCE 
provides an ideal tool for examining the precise behav- 
ior of the crossing point. We show its position, and the 
value of C v at the crossing, in Figs.^c) and^d), respec- 
tively. We find that, for fermions, the temperature and 
the specific heat value of the crossing point between two 
consecutive values of U is nearly independent of U for 
U < 8 (T cross - - 1.75, and C£ ross - - 0.36), with changes 
of roughly 5%. These variations are slightly larger for 
2HCBs in the same range of interactions. 

Away from half filling, where the ground state of nei- 
ther system is known, the trends in the deviation of finite- 
temperature properties of fermions and 2HCBs is differ- 
ent from the one reported so far. As an example, we 
show in Fig. [3] the entropy (for U — 4, 8, and 16) vs T 
for n = 0.85, and compare each curve with the one for 
n = 1.00. Interestingly, in the weak-coupling regime, the 
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FIG. 4. (Color online) Equation of state and compressibility 
(dn/dfi) vs /i of fermions and 2HCBs for (a) U = 9 and (b) 
[/ = 14 at T = 0.47. Thick lines show n and thin lines 
show the compressibility. At zero temperature, the density 
of fermions is pinned at unity for \fj,\ < A, where A is the 
charge gap, regardless of U. For bosons, the system has no 
gap for U <U C . Here, Wynn resummations with three cycles 
of improvement are used, and only the last order within the 
convergence region is shown. 



entropy does not change significantly at the accessible 
temperatures for either particle statistics as one dopes 
the system away from half filling. For U = 8 [Fig. EJJb)], 
the entropy for the fermionic system with 15% doping 
is larger than the half-filled value for 0.3 < T < 2, 
while it does not change nearly as much with doping for 
2HCBs for T > 0.9. By increasing the interaction to 
U = 16, the agreement between the entropy of fermions 
and 2HCBs away from half filling improves, yet the devi- 
ations between the two for n = 0.85 start at much higher 
temperatures in comparison to the half filled case [see 
Fig. (3Jc)]. These observations suggest that, even in the 
strong-coupling regime, the two systems away from half 
filling have fundamentally different phases. 



Equation of state and compressibility 

Further insight on the phases of the two systems at low 
temperatures can be gained by studying their equations 
of state and compressibilities. They are also of great 
interest to optical lattice experiments since those exper- 
iments are done in the presence of a confining potential 
that imposes a spatially varying chemical potential on 
the system. So, different regions in the trap correspond 
to different densities. In Fig. 21 we show n and dn/dfi 
vs fi at T — 0.47 for U — 9 and U = 14, which are be- 
low and above the critical interaction value for 2HCBs. 
From QMC calculations it is known that, for U < U c , the 
ground state of 2HCBs does not have a charge density 
gap [h|. One can also infer from Fig. [3] that, regardless 
of the interaction strength, 2HCBs always have a smaller 
Mott gap at zero temperature than the fermions. This, 
in turn, implies that around half filling, the compressibil- 
ity is always greater for bosons than for fermions, at any 
given temperature. 
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FIG. 5. (Color online) Double occupancy vs temperature 
for different values of U at half filling (top panels) and for 
n — 0.85 (bottom panels). In the weak-coupling regime with 
(7 = 4 [(a) and (d)], D for 2HCBs increases significantly by 
lowering the temperature below 1. Thick (color) lines are re- 
sults from the last order and thin (black) lines are the results 
from the one-to-last order of the NLCEs after resummations. 



Double occupancy 

A slight upturn in the double occupancy (D) as one de- 
creases T for the half-filled strongly-correlated fermionic 
system is a known phenomenon that has been attributed 
to the increase in virtual hoppings between allowed near- 
est neighbor sites due to the enhancement of short-range 
AF correlations [3]|. In fact, the onset of this increase, 
which can be measured in the experiments, may serve as 
a universal probe for large AF correlations |36( ■ It can be 
shown to lead to adiabatic cooling by increasing the inter- 
action strength [3(J , which is of great interest to optical 
lattice experiments. We have recently shown that such 
an increase in the double occupancy also occurs away 
from half-filling. Hence, in optical lattice experiments 
one also needs to make sure the density in most of the 
system is around half-filling in order for any increase in 
the double occupancy to be associated with the onset of 
antifcrromagnctism [2l| . 

We find interesting trends in the double occupancy of 
2HCBs when compared to fermions, especially for weak 
interactions. As shown in Fig. [5j and unlike in the 
fermionic case, the double occupancy of 2HCBs at and 
away from half filling increases sharply below T ~ 1 
for weak interactions, e.g., U = 4 in Figs. [Sfa) and 
EJd), while such a large difference between the results 
for the two particle statistics is absent for large U. The 
sharp low-T rise in D indicates that adiabatic cooling 
starting from the weakly interacting limit is efficient for 
2HCBs The double occupancy for 2HCBs even 

reaches the uncorrelated value of 1/4 in the half-filled 
case for U = 4 and T ~ 0.2, in stark contrast to D of the 
fermionic system. This is consistent with the absence of 
a Mott insulator for 2HCBs in that parameter region. 
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FIG. 6. (Color online) Nearest-neighbor spin correlations vs 
temperature at half filling for (a) (7 = 4, (b) (7 = 8, and 
(c) U — 16. (d)-(f) Uniform spin susceptibility at half filling 
vs temperature. Like for fermions, the spin susceptibility for 
2HCBs peaks at a characteristic temperature T*. The inset in 
(a) shows the exponential decrease of the difference between 
S zz of 2HCBs and fermions at fixed temperatures by increas- 
ing ?7 in the strong-coupling regime. Lines are the same as in 
Fig. [5] 



Spin correlations and uniform susceptibility 



A related trend is also seen in the results for the NN 
spin correlations, S zz = \§2$SfSf +s )\, where the sum 
runs over the four nearest neighbors of site i. As de- 
picted in Fig.EJa), S zz for 2HCBs peaks around T = 0.5 
when (7 = 4 before becoming vanishingly small at lower 
temperatures, which is again consistent with the super- 
fluid nature of the ground state in this interaction re- 
gion. On the other hand, it has been shown that for the 
fermionic case, S zz at half filling grows monotonically 
by decreasing the temperature [32(, which is consistent 
with its antiferromagnetically ordered low-T phase. It 
would be interesting to examine S zz at T = for 2HCBs 
and across the phase transition between the 2SF and the 
SCF phases, where S zz is presumably not small. Here, 
the difference between results for fermions and 2HCBs 
sets in around T = 1 for U = 4 and, like all other ther- 
modynamic quantities at half filling, becomes smaller as 
U is increased [see Figs. H2a)-(c)]. As shown in the in- 
set of Figs. EJa), at fixed temperatures, this difference 
becomes exponentially small with increasing U in the 
strong-coupling regime. 

Another thermodynamic quantity that highlights the 
difference between fermions and 2HCBs in the weak- 
coupling regime is the uniform spin susceptibility (x). 
As seen in Fig. EJd)-(f), we find that the deviation of \ 
for 2HCBs from that of fermions is large at low temper- 
ature for U = 4 and becomes smaller as U increases to 
8 and 16. More importantly, whereas a previous QMC 
study by Paiva et al. [2(| has shown that the peak lo- 
cation in the fermionic case changes non-monotonically 
by increasing U (following the variations of the AF cor- 



relations in the system), for 2HCBs, the peak tempera- 
ture decreases monotonically by increasing the interac- 
tion strength. This is because, unlike for fermions, the 
peak in x f° r 2HCBs in the weak-coupling regime does 
not signify moment ordering, but rather the disappear- 
ance of well-defined moments. This can be understood 
from the fact that NN spin correlations also decrease 
around the same temperature [see Fig. [5]Ja)]. 



Attractive interactions 

It is known for the fermionic Hubbard model that 
there exists the following unitary particle-hole transfor- 
mation [53 that takes the repulsive Hubbard model to 
the attractive one: 



fit = d 



f 



jl 



j'f > 



(8) 



where Rj is the displacement vector of site j. To sec 
the effect of the transformation more clearly, it is easier 
to rewrite the Hamiltonian of Eq. ([1]) in the so-called 
particle-hole symmetric form: 



(9) 



where fx' = fj, — ^. The transformation in Eq. ((8]) leaves 
the hopping term, as well as n^, invariant, but changes 
Tin to 1 — nti- As a result, the Hamiltonian in Eq. Q is 
transformed to: 

# = -tJ^(&4r+H.C.) (10) 



where 



dLdi, 



At half filling ((J = 0), the only 
change from Eq. © will be the sign of the interaction 
U. Therefore, the energy spectral properties of the half- 
filled Fcrmi-Hubbard model, e.g., its specific heat, en- 
tropy, etc, are invariant under U —>—[/. The nature of 
the ground state, however, is profoundly different in the 
repulsive and the attractive models since this transfor- 
mation maps charge correlations to spin correlations and 
vice versa J38| , which means that long-range AF order is 
mapped to a chargc-density-wave one. 

A similar unitary transformation maps the repulsive 
Hubbard model for 2HCBs to an attractive one: 



b ir = d ir 
hi = d 



jl' 



(11) 



Therefore, similar to the fermionic case, the Hamiltonian 
is invariant under the change of sign of the interaction at 
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FIG. 7. (Color online) Specific heat of the repulsive Hubbard 
model for 2HCBs with and without a magnetic field, h — 0.5 
and 0, for (a) U = 8 and (b) U = 16. Thick (thin) lines 
are the results from the last (one-to-last) order of the NLCEs 
after numerical resummation. 

half filling. The same argument presented above for the 
nature of the ground state of the fermionic model also 
applies to 2HCBs for U > U c . In the regime U < U Cl we 
expect the 2SF phase to be the ground even with attrac- 
tive interactions since, like the repulsive case, the system 
presumably gains more energy through the condensation 
of each specie than by minimizing the interaction energy. 

As is clear from Eqs. © and (fTT]) . the attractive (re- 
pulsive) Hamiltonian away from half filling (// ^ 0) is 
equivalent to the repulsive (attractive) one in the pres- 
ence of a magnetic field in the z direction, h, a role that 
is played by the chemical potential in the attractive (re- 
pulsive) case. In Fig. [71 we show how the specific heat 
of the repulsive Hubbard model for 2HCBs is modified 
in the presence of such a field, by plotting C v of the at- 
tractive model at // = and 0.5 (h = and 0.5 for the 
repulsive case) . The low temperature region is not acces- 
sible to us at small U and, in Fig. UJa) for U — 8, only 
a small deviation around the high-temperature peak can 
be seen when the magnetic field is introduced. The re- 
sults for U — 16 in Fig. [TJb) show a suppression of the 
specific heat at T < 1, which are consistent with the fact 
that spin degrees of freedom emerge only in the latter 
temperature region. Results for the fermionic case show 
qualitatively the same behavior as for 2HCBs for those 
two values of U. 



Adiabatic cooling: 
implications for optical lattice experiments 

Previously, we mapped out the isentropic paths of 
the 2D fermionic Hubbard model in the extended 
temperature- interaction space [2l[. In Fig. Ufa), we 
present a similar diagram for 2HCBs. As expected from 
the large negative slope in the low-T double occupancy 
of 2HCBs for small U vs T (Fig. [5]), adiabatic cooling by 
increasing U in the weak-coupling regime is much more 
efficient in comparison to the fermionic case. With the 
entropy per particle of 0.6, T reduces roughly by a factor 
of 4 as the interaction increases from 1 to 24. However, 
the underlying physics in different regions of the diagram 




FIG. 8. (Color online) (a) Isentropic dependence of tempera- 
ture on the interaction for different values of the entropy for 
the 2HCB-Hubbard model at half filling. Also shown in this 
panel is the characteristic temperature, T* , which is the loca- 
tion of the peak in the uniform spin susceptibility, and T sp , 
which is the location of the peak in the NN spin correlations. 
The exact values of the latter are at lower temperatures than 
what is accessible to us for U < 4 and U > 6. (b) Isothermic 
curves of S zz vs U. (c) Isentropic curves of S zz vs U for the 
entropy per particle from 0.3 to ln(2). 



and the change in the strength of AF correlations is unlike 
that of fermions. As mentioned before, in the fermionic 
system, regardless of the value of the interaction, AF 
correlations are always enhanced by lowering the tem- 
perature. So, based on their dependence on U and how 
rapidly the temperature falls in the adiabatic cooling pro- 
cess, one could drive the fermionic system into the region 
with exponentially large AF correlations. The latter re- 
gion could be identified by the onset of a downturn in the 
uniform susceptibility (T*) |2ll ]. 

On the contrary, for 2HCBs, the peak in the uniform 
susceptibility does not signify large AF correlations for 
weak interactions [see Figs. Ufa) and[5Jd)]. We empha- 
size this by plotting in the same figure the location of 
the peak in S zz (T sp ) for a few values of U. Beyond 
U c , where the Mott insulating ground state not only has 
long-range ccy-ferromagnetic order, but also very large 
z-antiferromagnetic correlations, the maximum of S zz is 
likely at zero temperature. 

Results in Fig. [8jb) and Fig. [HJc) provide further in- 
sight into the behavior of S zz for different interaction 
strengths. By comparing these figures to their fermionic 
counterparts [Figs. 3(c) and 3(d) in Ref. [2l[, one can see 
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that unlike in the weak-coupling regime, where the results 
are qualitatively different for the two particle statistics, 
in the strong-coupling regime they are very similar [as 
shown in Fig. EJa) , the difference decreases exponentially 
with U] . As seen in Fig. EJb) , at fixed (and low) temper- 
atures S zz for 2HCBs sharply drops to small values by 
decreasing U. In Fig. [SJc), the isentropic curves of S zz 
are almost on top of each other for U < 7 for the en- 
tropics shown, and are expected to become vanishingly 
small at lower entropies that are not accessible to us. 
This is, again, unlike the trends in the fermionic model 
in the weak-coupling regime. The opposite is true, how- 
ever, for U > U c where, at the entropies accessible to us, 
S zz values of 2HCBs agree very well with their fermionic 
counterparts. 



So far, much lower temperatures have been achieved 
with bosons in optical lattices than with fermions. Em- 
ploying novel cooling techniques, researchers have been 
able to access temperatures as low as 1 nK with bosons 
deep in the Mott insulating regime [39[ . It has also been 
shown in experiments with bosonic mixtures and that the 
inter- and intra-specie interactions can be tuned using 
Feshbach resonance and other techniques @, H, OH • Our 
results show that thermodynamic properties and short- 
range AF correlations of 2HCBs are very (exponentially) 
close to those of the fermions for strong interactions and 
up to intermediate to low temperatures. Therefore, by 
engineering the strong interaction regime studied here, 
optical lattice experiments with two-component bosons 
could also provide a tool for simulating intermediate 
to low temperature correlations in the Fermi-Hubbard 
model. Further studies need to be done to explore the 
properties of 2HCBs away from half filling and examine 
their relevance to the fermionic case. 



IV. CONCLUSIONS 

We have examined the particle statistics dependence of 
the thermodynamic properties of the 2D Hubbard model 
by means of an exact method (NLCEs) that yields prop- 
erties such as entropy, specific heat, double occupancy, 
and spin correlations in the thermodynamic limit. The 
results are valid above a certain temperature that de- 
pends on the interaction and the filling factor. We consid- 
ered two-component fermions and 2HCBs and compared 
their properties at various temperatures and interaction 
strengths, up to three times the bandwidth. 

We have shown that for weak interactions, in the 
regime where the ground state of the two models at half 
filling have fundamentally different natures, the results 
for observables at finite temperature differ significantly 
for the two particle statistics starting at relatively high 
temperatures. In contrast, in the strong-coupling regime 
(beyond a critical interaction that separates a superfluid 
from a Mott insulator phase for 2HCBs), the agreement 
between the thermodynamic quantities of the two sys- 
tems, including short-range AF correlations, extends to 
much lower temperatures. We find that the trends in spin 
correlations in the z direction are similar for both particle 
statistics and that the results differ only by exponentially 
small values for strong interactions. This provides an ad- 
ditional tool to probe correlations of fermionic systems 
in optical lattice experiments by emulating them using 
two-specie bosons, which can generally be cooled down 
to lower temperatures. 
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